NASA  Contractor  Report  194992 
ICASE  Report  No.  94-81 


ICASE 

AN  ANALYSIS  OF  SPECTRAL  ENVELOPE-REDUCTION 
VIA  QUADRATIC  ASSIGNMENT  PROBLEMS 


Alan  George 
Alex  Pothen 


Contract  NASl-19480 
October  1994 


19950104  064 


Institute  for  Computer  Applications  in  Science  and  Engineering 
NASA  Langley  Research  Center 
Hampton,  VA  23681-0001 


operated  by  Universities  Space  Research  Association 

DISTmBUTXOi"STATaS^^^ 

Approved  for  public  release; 

Distribution  Unlimited 


AN  ANALYSIS  OF  SPECTRAL  ENVELOPE-REDUCTION 
VIA  QUADRATIC  ASSIGNMENT  PROBLEMS  * 

ALAN  GEORGEt  AND  ALEX  POTHENi 


Abstract.  A  new  spectral  algorithm  for  reordering  a  sparse  symmetric  matrix  to  reduce  its  envelope 
size  was  described  in  [2].  The  ordering  is  computed  by  associating  a  Laplacian  matrix  with  the  given  matrix 
and  then  sorting  the  components  of  a  specified  eigenvector  of  the  Laplacian.  In  this  paper  we  provide  an 
analysis  of  the  spectral  envelope  reduction  algorithm.  We  describe  related  1-  and  2-sum  problems;  the 
former  is  related  to  the  envelope  size,  while  the  latter  is  related  to  an  upper  bound  on  the  work  involved  in 
an  envelope  Cholesky  factorization  scheme.  We  formulate  the  latter  two  problems  as  quadratic  assignment 
problems,  and  then  study  the  2-sum  problem  in  more  detail.  We  obtain  lower  bounds  on  the  2-sum  by 
considering  a  projected  quadratic  assignment  problem,  and  then  show  that  finding  a  permutation  matrix 
closest  to  an  orthogonal  matrix  attaining  one  of  the  lower  bounds  justifies  the  spectral  envelope  reduction 
algorithm.  The  lower  bound  on  the  2-sum  is  seen  to  be  tight  for  reasonably  “uniform”  finite  element  meshes. 
We  also  obtain  asymptotically  tight  lower  bounds  for  the  envelope  size  for  certain  classes  of  meshes. 
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1.  Introduction.  A  novel  spectral  algorithm  to  reduce  the  envelope  of  a  sparse,  sym¬ 
metric  matrix  was  described  in  a  companion  paper  [2].  The  algorithm  associates  a  discrete 
Laplacian  matrix  with  the  given  symmetric  matrix,  and  then  computes  a  reordering  of  the 
matrix  by  sorting  the  components  of  an  eigenvector  corresponding  to  the  smallest  nonzero 
Laplacian  eigenvalue.  The  results  in  [2]  show  that  the  spectral  algorithm  can  obtain  sig¬ 
nificantly  smaller  envelope  sizes  compared  to  other  currently  used  algorithms.  All  previous 
envelope-reduction  algorithms  (known  to  us),  e.g.,  the  reverse  Cuthill-McKee  (RCM)  algo¬ 
rithm,  the  Gibbs-Poole-Stockmeyer  (GPS)  algorithm,  the  Gibbs-King  (GK)  algorithm,  and 
the  Sloan  algorithm  [3,  14,  15,  23,  36],  are  combinatorial  in  nature,  employing  some  variant 
of  breadth-first-search  to  compute  the  ordering.  In  contrast,  the  spectral  algorithm  is  an 
algebraic  algorithm  whose  good  envelope-reduction  properties  are  somewhat  intriguing  and 
poorly  understood. 

In  this  paper  we  attempt  to  provide  a  raison  d’etre  for  the  spectral  envelope- reduction 
algorithm.  We  describe  problems  related  to  envelope-reduction  called  the  1-  and  2-sum  prob¬ 
lems,  and  then  formulate  these  latter  problems  as  quadratic  assignment  problems  (QAPs). 
We  show  that  the  QAP  formulation  of  the  2-sum  enables  us  to  obtain  lower  bounds  on  the 
2-sum  (and  related  envelope  parameters)  based  on  the  Laplacian  eigenvalues.  The  lower 
bounds  seem  to  be  quite  tight  for  finite  element  problems  whose  mesh  points  are  nearly  all 
of  the  same  degree.  Further,  we  show  that  a  closest  permutation  matrix  to  an  orthogonal 
matrix  that  attains  the  lower  bound  is  obtained  by  sorting  the  second  Laplacian  eigenvec¬ 
tor  components  in  monotonically  increasing  or  decreasing  order.  This  provides  stronger 
justification  for  the  spectral  envelope-reducing  algorithm  than  has  been  provided  earlier. 

The  computational  results  in  [2]  indicate  that  the  spectral  envelope  reduction  algorithm 
is  even  more  effective  in  reducing  the  work  in  an  envelope  factorization  algorithm,  due  to  the 
observed  quadratic  dependence  of  the  factorization  time  on  the  envelope  size.  The  analysis 
in  this  paper,  in  terms  of  the  2-sum  problem,  explains  this  phenomenon  as  well. 

Although  initially  envelope-reducing  orderings  were  developed  for  use  in  envelope  schemes 
for  sparse  matrix  factorization,  these  orderings  have  been  used  in  the  past  few  years  in  sev¬ 
eral  other  applications.  The  RCM  ordering  has  been  found  to  be  an  effective  preordering 
in  computing  incomplete  factorization  preconditioners  for  preconditioned  conjugate- gradient 
methods  [4,  6].  Envelope-reducing  orderings  have  been  used  in  frontal  methods  for  sparse 
matrix  factorization  [7].  Such  orderings  have  also  been  used  in  paralleFInatrix- vector  multi¬ 
plication  and  tridiagonalization  of  sparse  symmetric  matrices. 

The  wider  applicability  of  envelope-reducing  orderings  prompts  us  to  take  a  fresh  look 
at  the  reordering  algorithms  currently  available,  and  to  develop  new  ordering  algorithms. 
Spectral  envelope- reduction  algorithms  seem  to  be  attractive  in  this  context,  since  they 

(i)  compare  favorably  with  existing  algorithms  in  terms  of  the  quality  of  the  orderings  [2], 

(ii)  extend  easily  to  problems  with  weights,  e.g.,  finite  element  meshes  arising  from  dis¬ 
cretizations  of  anisotropic  problems,  and 

(iii)  are  fairly  easily  parallelizable. 

Spectral  algorithms  are  usually,  but  not  always,  more  expensive  than  the  other  algorithms 
currently  available,  but  since  the  envelope-reduction  problem  requires  only  one  eigenvector 
computation  (to  a  few  decimal  digits  of  precision),  we  believe  the  costs  are  not  impractically 
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high.  Improvements  to  reduce  the  costs  are  being  looked  into  as  well.  We  focus  primarily  on 
the  class  of  finite  element  meshes  arising  from  discretizations  of  partial  differential  equations. 
Our  goals  in  this  project  are  to  develop  efficient  software  implementing  our  algorithms,  and 
to  prove  results  about  the  quality  of  the  orderings  generated. 

The  projection  approach  for  obtaining  lower  bounds  of  a  QAP  is  due  to  Hadley,  Rendl, 
and  Wolkowicz  [17],  and  this  approach  has  been  applied  to  the  graph  partitioning  problem 
by  the  latter  two  authors  [34].  In  earlier  work  a  spectral  approach  for  the  graph  (matrix) 
partitioning  problem  has  been  employed  to  compute  a  spectral  nested  dissection  ordering 
for  sparse  matrix  factorization,  for  partitioning  computations  on  finite  element  meshes  on  a 
distributed-memory  multiprocessor  [19,  31,  32,  35],  and  for  load-balancing  parallel  compu¬ 
tations  [20].  The  spectral  approach  has  also  been  used  to  find  a  pseudo-peripheral  node  [16]. 
Juvan  and  Mohar  [21,  22]  have  provided  a  theoretical  study  of  the  spectral  algorithm  for 
reducing  p-sums,  where  p  =  1,  2,  and  oo,  and  Helmberg  et  al.  [18]  obtain  spectral  lower 
bounds  on  the  bandwidth.  A  survey  of  some  of  these  earlier  results  may  be  found  in  [29]. 
Paulino  et  al.  [30]  have  also  considered  the  use  of  spectral  envelope-reduction  for  finite  ele¬ 
ment  problems. 

The  following  is  an  outline  of  the  rest  of  this  paper.  In  Section  2  we  describe  various 
parameters  of  a  matrix  associated  with  its  envelope,  introduce  the  envelope  size  and  envelope 
work  minimization  problems,  and  the  related  1-  and  2-sum  problems.  Let  A  denote  the 
maximum  number  of  offdiagonal  nonzeros  in  a  row  or  column  of  the  given  matrix  (the 
maximum  degree  of  a  vertex  in  the  adjacency  graph  of  the  matrix).  We  prove  that  the 
minimum  1-sum  is  bounded  below  by  the  minimum  envelope  size,  and  bounded  above  by 
A  times  the  minimum  envelope  size.  A  similar  result  holds  for  the  minimum  2-sum  and  a 
bound  on  the  work  in  an  envelope  Cholesky  factorization.  We  compute  upper  and  lower 
bounds  for  the  envelope  of  a  sparse  symmetric  matrix  in  terms  of  the  eigenvalues  of  the 
Laplacian  matrix  in  Section  3.  The  popular  RCM  ordering  is  obtained  by  reversing  the 
Cuthill-McKee  (CM)  ordering,  and  usually  the  RCM  ordering  has  smaller  envelope  size 
and  work  than  the  CM  ordering.  We  prove  that  reversing  an  ordering  can  improve  or 
impair  the  envelope  size  by  at  most  a  factor  A,  and  the  envelope  work  by  at  most  A^.  In 
Section  4,  we  formulate  the  2-  and  1-sum  problems  as  quadratic  assignment  problems.  We 
obtain  lower  and  upper  bounds  for  the  2-sum  problem  in  terms  of  the  eigenvalues  of  the 
Laplacian  matrix  in  Section  5  by  means  of  a  projection  approach  that  relaxes  a  permutation 
matrix  to  an  orthogonal  matrix  with  row'  and  column  sums  equal  to  one.  We  also  show 
how  the  lower  bound  may  be  strengthened  by  permitting  diagonal  perturbations  to  the 
Laplacian  matrix.  We  justify  the  spectral  envelope-reduction  algorithm  in  Section  6  by 
proving  that  a  closest  permutation  matrix  to  an  orthogonal  matrix  attaining  the  lower  bound 
for  the  2-sum  is  obtained  by  permuting  the  second  Laplacian  eigenvector  in  monotonically 
increasing  or  decreasing  order.  In  Section  7  we  employ  the  lower  bounds  we  have  obtained 
to  study  the  asymptotic  growth  of  envelope  size  and  an  estimate  of  the  envelope  work  for 
2-  and  3-  dimensional  finite  element  meshes  of  bounded  degree.  We  show  that  graphs  with 
good  separators  have  small  envelope  parameters  as  well,  by  considering  a  modified  nested 
dissection  ordering.  We  present  computational  results  in  Section  8  to  illustrate  that  the  2- 
sums  obtained  by  the  spectral  reordering  algorithm  can  be  close  to  optimal  for  finite  element 
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meshes  whose"  nodes  have  almost  equal  degrees.  We  make  some  concluding  remarks  in 
Section  9.  The  Appendix  contains  some  lower  bounds  for  the  more  general  p-sum  problem, 
where  1  <  p  <  oo. 

2.  A  menagerie  of  envelope  problems. 

2.1.  The  envelope  of  a  matrix.  Let  A  be  an  n  x  n  symmetric  matrix  with  elements 
Cij,  whose  diagonal  elements  are  nonzero.  Various  parameters  of  the  matrix  A  associated 
with  its  envelope  are  defined  below. 

We  denote  the  column  indices  of  the  nonzeros  in  the  lower  triangular  part  of  the  zth  row 
by 


row{i)  =  {j  :  aij  ^  0  and  1  <j  <  *}• 

For  the  ith  row  of  A  we  define 

fi{A)  =  min{j  :  j  e  ro^o(^)},  and 
ri{A)  =  i  -  fi{A). 

Here  fi{A)  is  the  column  index  of  the  first  nonzero  in  the  Ah  row  of  A  (by  our  assumption 
of  nonzero  diagonals,  1  <  fi  <  i),  and  the  parameter  ri(A)  is  the  row-width  of  the  Ah  row 
of  A.  The  bandwidth  of  A  is  the  maximum  row-width 


bw{A)  =  max{r,(A)  :  i  =  1, . . . ,  n}. 

The  envelope  of  A  is  the  set  of  column  indices 

Env{A)  =  :  fi{A)  <j<i,  z  =  l,...,n}. 

For  each  row,  these  indices  lie  in  an  interval  beginning  with  the  column  index  of  the  first 
nonzero  element  and  ending  with  (but  not  including)  the  index  of  the  diagonal  nonzero 
element. 

We  denote  the  size  of  the  envelope  by  Esize{A)  =  \Env{A)\.  (This  is  also  called  the 
profile  of  A.)  The  work  in  the  Cholesky  factorization  of  A  that  employs  an  envelope  storage 
scheme  can  be  bounded  by  (1/2)  X)"=2  +  3).  This  bound  can  differ  significantly  from 

the  actual  work  required  on  some  problems,  but  for  nearly  uniform  finite  element  problems 
without  appendages,  it  is  usually  a  good  measure  of  the  work.  (The  computational  results 
on  the  work  in  an  envelope  factorization  scheme  in  [2]  indicate  this  is  a  good  approximation 
as  well.)  Hence  hereafter  we  will  use 


Wbound{A)  —  rf 

i=l 

as  a  measure  of  the  work  in  such  a  factorization. 

A  3  X  3  7-point  grid  and  the  nonzero  structure  of  the  corresponding  matrix  A  are  shown 
in  Figure  2.1.  A  ‘  ’  indicates  a  nonzero  element,  and  a  ‘  ’  indicates  a  zero  element 

that  belongs  to  the  lower  triangle  of  the  envelope  in  the  matrix.  The  row-widths  given  in 
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Fig.  2.1.  An  ordering  of  7-poini  grid  and  the  corresponding  matrix, 
is  indicated  by  marking  zeros  within  it  by  asterisks. 


The  lower  triangle  of  the  envelope 
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Table  2.1 


Row-widths  and  column-widths  of  the  matrix  in  Figure  1. 


Table  2.1  are  easily  verified  from  the  structure  of  the  matrix.  The  envelope  size  is  obtained 
by  summing  the  row-widths,  and  is  equal  to  24.  (Column-widths  c,-  are  defined  later  in  this 
section.) 

The  values  of  these  parameters  strongly  depend  on  the  choice  of  an  ordering  of  the  rows 
and  columns.  Hence  we  consider  how  these  parameters  vary  over  symmetrically  permuted 
matrices  P^AP,  where  P  is  a  permutation  matrix.  We  define  Esizem.in{7^),  the  minimum 
envelope  size  of  A,  to  be  the  minimum  envelope  size  among  all  permuted  matrices  P^ AP. 
The  quantities  Wboundmini^)  and  hwmin{A)  are  defined  in  similar  fashion.  Minimizing  the 
envelope  size  and  the  bandwidth  of  a  matrix  are  NP-complete  problems  [25],  and  minimizing 
the  work  bound  is  likely  to  be  intractable  as  well.  So  one  has  to  settle  for  heuristic  orderings 
to  reduce  these  quantities. 

It  is  helpful  to  consider  a  ‘column-oriented’  expression  for  the  envelope  size  for  obtaining 
a  lower  bound  on  this  quantity  in  Section  3.  The  width  of  a  column  j  of  A  is  the  set  of  all 
row  indices  in  the  jth  column  of  the  envelope  of  A.  In  other  words, 

Cj{A)  =  \{k  :  k  >  j, and  <  j^akt  0}|. 
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(This  is  also  called  the  jth  front-width.)  It  is  then  easily  seen  that  the  envelope  size  is 

n 

(2.1)  Esize{A)  = '^Cj. 

j=i 

The  work  in  an  envelope  factorization  scheme  is  given  by 

n 

Ework(A)  =  (1/2)  ^2  +  3). 

Hereafter,  we  will  ignore  the  linear  term  in  Cj  in  computing  the  envelope  work.  The  column- 
widths  of  the  matrix  in  Figure  2.1  are  given  in  Table  2.1.  These  concepts  and  their  inter¬ 
relationships  were  considered  earlier  by  Liu  and  Sherman  [27],  and  are  also  discussed  in  the 
text  books  [5,  13]. 

The  envelope  parameters  can  also  be  defined  with  respect  to  the  adjacency  graph  G  — 
{V,E)  of  A.  Denote  nbr{v)  —  {u}  U  adj{v).  In  terms  of  the  graph  G  and  an  ordering  a  of 
its  vertices,  we  can  define 

r(u,  o;)  =  max{a(u)  —  a(w)  :  w  G  nbr{v),  a(w)  <  Q!(u)}, 

P(G,a)  =  ma,x{a{v)  —  a(w)  :  {v,w)  ^  E}. 

Hence  we  can  write  the  envelope  size  and  work  associated  with  an  ordering  a  as 

Esize{G,a)  =  ^  r(u)  =  max{a(u)  —  a(w)  :  w  G  nbr{v).,a{w)  <  q(u)} 

v^V  v^V 

Wbound{G,a)  =  '^r{v)^  —  'Y^ma,x{{a{v)  —  a{w))‘^  :  w  €.  nbr{v),a{w)  <  a{v)}. 

v^V  v^V 

The  goal  is  to  choose  a  vertex  ordering  a  :  V  i-j-  {l,...,n}  to  minimize  one  of  the  param¬ 
eters  described  above.  We  denote  by  Esizemin{G)  {Wboundmin{G))  the  minimum  value  of 
Esize{G,  a)  (Wbound[G,  a))  over  all  orderings  a.  The  reader  can  compute  the  envelope  size 
of  the  numbered  graph  in  Figure  2.1,  using  the  definition  given  in  this  paragraph,  to  verify 
that  Esize{G)  =  24. 

The  jth  front- width  has  an  especially  nice  interpretation  if  we  consider  the  adjacency 
graph  G  =  (V,  E)  of  A.  Let  the  vertex  corresponding  to  a  column  j  of  A  be  numbered  Vj  so 
that  V  =  {ui, . . . ,  u„},  and  define  Vj  =  {ui, . .  Denote  adj{X)  =  {Vy^xO'dj{v))  \  X,  for 

a  subset  of  vertices  X.  Then  Cj{A)  —  \adj{yj)\. 

To  illustrate  the  dependence  of  the  envelope  size  on  the  ordering,  we  include  in  Figure  2.2 
an  ordering  that  leads  to  a  smaller  envelope  size  for  the  7-point  grid.  Again,  a  ‘  ’  indicates 
a  nonzero  element,  and  a  ‘  ’  indicates  a  zero  element  that  belongs  to  the  lower  triangle  of 

the  envelope  in  the  matrix.  This  ordering  by  ‘diagonals’  yields  the  optimal  envelope  size  for 
the  7-point  grid  [24]. 

2.2.  1-  and  2-sum  problems.  It  will  be  helpful  to  consider  quantities  related  to  the 
envelope  size  and  envelope  work,  the  1-sum  and  the  2-sum.  We  write  the  envelope  size  and 
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10' - ^ ^ ^ - ' 

0  2  4  6  8  10 

9  Esize=19 

Fig.  2.2.  Another  ordering  of  a  7-point  grid  and  the  corresponding  matrix.  Again  the  lower  triangle  of 
the  envelope  is  indicated  by  marking  the  zeros  within  it  by  asterisks. 

1-sum,  and  the  envelope  work  and  the  2-sum  in  a  way  that  shows  their  relationships: 


n 


(2.2) 

Esize(A)  =  ^  max  (i  —  j), 
i=i  jerow{i) 

(2.3) 

^i(^)  =  E  E 

j€row{i) 

(2.4) 

n 

Whound{A)  =  ^  max  {i  —  jY, 
i=i  jerow{i) 

(2.5) 

=  E  E 

jerow{i) 

The  parameters  Ui^miniA)  and  are  the  minimum  values  of  these  parameters  over 

all  permuted  matrices  AP. 

More  generally,  for  real  1  <  p  <  oo,  we  define  the  p-sum  to  be 

jerow{i) 

The  1-sum  problem  (p 
and  the  limiting  case  p 

=  1)  has  also  been  called  the  optimal  linear  arrangement  problem, 
=  oo  corresponds  to  the  bandwidth  problem;  both  these  problems  are 

NP-complete. 

We  now  consider  the  relationships  between  the  envelope  size  problem  and  the  1-sum 
problem,  and  between  the  envelope  work  problem  and  the  2-sum  problem.  Let  A  denote 
the  maximum  number  of  offdiagonal  nonzeros  in  a  row  of  A.  (This  is  the  maximum  vertex 
degree  in  the  adjacency  graph  of  ^4.) 
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Theorem  2.1.  The  minimum  values  of  the  envelope  size,  envelope' work  in  the  Cholesky 
factorization,  1-sum,  and  2-sum  of  a  symmetric  matrix  A  are  related  by  the  following  in¬ 
equalities: 


(2.6) 

< 

<^l,mjn(A)  ^  AEsiZCiYiini^Ad^', 

(2.7) 

Whoundmin{A) 

< 

^2,mm(^)  ^  A'Wb0Undmin{A)-, 

(2.8) 

< 

^  “sj  E 

Proof  We  begin  by  proving  (2.7).  Our  strategy  will  be  to  first  prove  the  inequalities 

Wbound{A)  <  (tKA)  <  AWbound{A), 


and  then  to  obtain  the  required  result  by  considering  two  different  permutations  of  A. 

The  bound  Wbound{A)  <  crl(A)  is  immediate  from  equations  (2.4)  and  (2.5).  If  the 
inner  sum  in  the  latter  equation  is  bounded  from  above  by 

A  max  (i  —  jY, 

j£row{i) 


then  we  get  AWbound{A)  as  an  upper  bound  on  the  2-sum. 

Now  let  Xi  be  a  permutation  matrix  such  that  Ai  =  AXi,  and  'Wbound{Ai)  = 
'Wboundmini'^)-  Then  we  have 

<^2,mm(^)  ^  <^2(^1)  ^  AWbound{Ai)  =  AWboundmin{A). 

Further,  let  X2  be  a  permutation  matrix  such  that  A2  =  A’^AA’2,  and  cr|(A2)  =  c>'2^mm(^)- 
Again,  we  have 

Wboundmin{A)  <  Wbound{A2)  <  C»’2(^2)  =  (^2,min{^)‘ 


We  obtain  the  result  by  putting  the  last  two  inequalities  together. 

We  omit  the  proof  of  (2.6)  since  it  can  be  obtained  by  a  similar  argument,  and  proceed 
to  prove  (2.8).  The  first  inequality  cr2(A)  <  cfx{A)  holds  since  the  p-norm  of  any  real  vector 
is  a  decreasing  function  of  p.  The  second  inequality  is  also  standardj^_ since  it  bounds  the 
1-norm  of  a  vector  by  means  of  its  2-norm.  This  result  was  obtained  earlier  by  Juvan  and 
Mohar  [22];  we  include  its  proof  for  completeness.  Applying  the  Cauchy-Schwarz  inequality 
to  crj(A)  we  have 


V 

jerow{i) 


A  \ 

(  n  \ 

<  j 

E  E  1 

E  E  (*-7? 

jerow{i)  ) 

ierow{i)  / 

\E\al{A). 


We  obtain  the  result  by  considering  two  orderings  that  achieve  the  minimum  1-  and  2-sums. 

□ 
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3.  Bounds  for  envelope  size.  In  this  section  we  present  some  bounds  for  the  envelope 
size  problem.  We  will  require  some  background  on  the  Laplacian  matrix. 

3.1.  The  Laplacian  matrix.  The  Laplacian  matrix  Q{G)  of  a  graph  G  is  the  n  x  n 
matrix  D  —  M,  where  D  is  the  diagonal  degree  matrix  and  M  is  the  adjacency  matrix  of  G. 
If  G  is  the  adjacency  graph  of  a  symmetric  matrix  A,  then  we  could  define  the  Laplacian 
matrix  Q  directly  from  A: 


'  -1 

if  i  /  j  and  Oij 

gij  —  ■' 

0 

if  i  ^  j  and  aij 

Y^]=i  \gij\ 

if  i  =  j. 

Note  that 

^Qx  = 

x^Dx  —  gy^Mx 

(3.1) 

= 

Y,  -  xjf- 

j<> 

a,j^0 

The  eigenvalues  of  Q(G)  are  the  Laplacian  eigenvalues  of  G,  and  we  list  them  as  Ai((5)  < 
A2(Q)  <  •  •  •  <  An((5).  An  eigenvector  corresponding  to  Xk(Q)  will  be  denoted  by  x*.,  and  will 
be  called  a  kth.  eigenvector  of  Q.  It  is  well-known  that  Q  is  a  singular  M-matrix,  and  hence 
its  eigenvalues  are  nonnegative.  Thus  Ai((5)  =  0,  and  the  corresponding  eigenvector  is  any 
nonzero  constant  vector  c.  If  G  is  connected,  then  Q  is  irreducible,  and  then  X2{Q)  >  0;  the 
smallest  nonzero  eigenvalues  and  the  corresponding  eigenvectors  have  important  properties 
that  make  them  useful  in  the  solution  of  various  partitioning  and  ordering  problems.  These 
properties  were  first  investigated  by  Fiedler  [8,  9];  as  discussed  in  Section  1,  more  recently 
several  authors  have  studied  their  application  to  such  problems. 

3.2.  Laplacian  bounds  for  the  envelope  size.  It  would  be  helpful  to  work  with 
the  ‘column-oriented’  definition  of  the  envelope  size.  Let  the  vertex  corresponding  to  a 
column  j  of  A  be  numbered  vj  in  the  adjacency  graph  so  that  V  =  {ui,...,^^},  and  let 
Vj  =  {ui,.  ..,Uj}.  Recall  that  the  width  of  a  column  j  of  A  is  Cj{A)  =  \adj{Vj)\^  and  that 
the  envelope  size  is 


n 

Esize{A)  —  c,. 

j=i 


Recall  also  that  A  denotes  the  maximum  degree  of  a  vertex,  and  that  given  a  set  of  vertices 
S',  we  denote  by  S{S)  the  set  of  edges  with  one  endpoint  in  S  and  the  other  in  R  \  S'. 

We  make  use  of  the  following  elementary  result,  where  the  lower  bound  is  due  to  Alon 
and  Milman  [1]  and  the  upper  bound  is  due  to  Juvan  and  Mohar  [22]. 

Lemma  3.1.  Let  S  CV  be  a  subset  of  the  vertices  of  a  graph  G.  Then 


X2iQ) 


|5||y\5| 


n 


<  |^(5)|  <  A„(Q) 


|5||y\5| 


n 


□ 
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Theorem  3.2.  The  envelope  size  of  a  symmetric  matrix  A  can  be  bounded  in  terms  of 
the  eigenvalues  of  the  associated  Laplacian  matrix  as 

-  1)  <  Esi.e(A)  <  -  1). 

Proof  From  Lemma  3.1, 

n 

Now  Cj{A)  =  \adj{Vj)\  >  |^(Vj)|/A;  substituting  the  lower  bound  for  |^(Vj)|,  and  summing 
this  latter  expression  over  all  we  obtain  the  lower  bound  on  the  envelope  size. 

The  upper  bound  is  obtained  by  using  the  inequality  Cj{A)  <  |^(Vj)|  with  the  upper 
bound  in  Lemma  3.1.  □ 

Cuthill  and  McKee  [3]  proposed  one  of  the  earliest  ordering  algorithms  for  reducing  the 
envelope  size  of  a  sparse  matrix.  George  [12]  discovered  that  reversing  this  ordering  often 
leads  to  a  significant  reduction  in  envelope  size  and  work.  Since  then  the  reverse- Cuthill- 
McKee  (RCM)  ordering  has  become  one  of  the  most  popular  envelope  size  reducing  orderings. 
However,  we  do  not  know  of  any  published  quantitative  results  on  the  improvement  that  may 
be  expected  by  reversing  an  ordering,  and  here  we  present  the  first  such  result.  For  degree- 
bounded  finite  element  meshes,  no  asymptotic  improvement  is  possible;  the  parameters  are 
improved  only  by  a  constant  factor. 

Theorem  3.3.  Reversing  the  ordering  of  a  sparse  symmetric  matrix  A  can  change 
(improve  or  impair)  the  envelope  size  by  at  most  a  factor  A,  and  the  envelope  work  by  at 
most  where  A  is  the  maximum  number  of  off  diagonal  nonzeros  in  a  row  (column)  of  A. 

Proof  Let  Vj  denote  the  vertex  in  the  adjacency  graph  corresponding  to  the  jfth  column 
of  A  (in  the  original  ordering)  so  that  the  jfth  column  width  Cj{A)  =  |ady(Lj)|,  where  Vj  = 
{uj, . . . ,  Uj}.  Let  A  denote  the  permuted  matrix  obtained  by  reversing  the  column  and  row 
ordering  of  A.  We  have  the  inequality 

Ci(A)  =  \aij(V,)\  <  |«(yi)|  <  A\adj(V  \  y,)|  =  Ac,.i(A). 

Since  Esize{A)  =  J2j=i  Cj(A),  summing  this  inequality  over  j  from  one  to  n,  we  obtain 
Esize{A)  <  AEsize{A).  By  symmetry,  the  inequality  Esize{A)  <~AEsize{A)  holds  as 
well. 

The  inequality  on  the  envelope  work  follows  by  a  similar  argument  from  the  equation 

Ework{A)  =  ° 

4.  Quadratic  assignment  formulation  of  2-  and  1-sum  problems.  We  formulate 
the  2-  and  1-sum  problems  as  quadratic  assignment  problems  in  this  section. 

4.1.  The  2-sum  problem.  Let  the  vector  p  —  I  2  n  and  let  o  be  a 
permutation  vector,  i.e.,  a  vector  whose  components  form  a  permutation  of  1,  . . .,  n.  We 
may  write  a  =  Xp,  where  A  is  a  permutation  matrix  with  elements 


I  1,  if  j  =  a{i) 
I  0,  otherwise 
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It  is  easily  verified  that  the  a{j))  element  of  the  permuted  matrix  X'^  AX  is  the  element 

aij  of  the  unpermuted  matrix  A.  Let  B  =  P  then  =  ij.  We  denote  the  set  of  all 
permutation  vectors  with  n  components  by 

We  write  the  2-sum  as  a  quadratic  form  involving  the  Laplacian  matrix  Q. 


=  Bf  Si  (“(')  “ 

—  "  a(j)<a(0 

=  min 

0€S„ 


~  " i=l j=l 


The  transformation  from  the  second  to  the  third  line  makes  use  of  (3.1). 

This  quadratic  form  can  be  expressed  as  a  quadratic  assignment  problem  by  substituting 
ba(i),c,(j)  =  o(?Xi): 

n  n 

minoXa  =  mm  ^  gy 

—  "  —  "  «=i  j=i 

There  is  also  a  trace  formulation  of  the  QAP  in  which  the  variables  are  the  elements  of 
the  permutation  matrix  X.  We  obtain  this  formulation  by  substituting  Xp  for  a.  Thus 

min  cx^Qa  =  mmp'^ X^ QXp. 

We  may  consider  the  last  scalar  expression  as  the  trace  of  a  1  x  1  matrix,  and  then  use  the 
identity  tr  MN  =  tr  NM  to  rewrite  the  right-hand-side  of  the  last  displayed  equation  as 

mintr  QXpp^X^  =  mintr  QXBX^ . 

X  X 

This  is  a  quadratic  assignment  problem  since  it  is  a  quadratic  in  the  unknowns  which 
are  the  elements  of  the  permutation  matrix  X.  The  fact  that  .6  is  a  rank-one  matrix  leads 
to  great  simplifications  and  savings  in  the  computation  of  good  lower  bounds  for  the  2-sum 
problem. 

4.2.  The  1-sum  problem.  Let  M  be  the  adjacency  matrix  of  a  given  symmetric  matrix 
A  and  let  S  denote  a  ‘distance  matrix’  with  elements  Sij  =  |i  —  j\,  both  of  order  n.  Then 

o'i,min{A)  =  mmai{X'^AX) 

X 

—  "  «(,)<«(.) 

n  n 

=  (1/2)™?  EE  *5a(z),or(j) 

i=l  j=l 

-  (l/2)minfr  A/X5X. 
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Unlike  the  2-sum,  the  matrices  involved  in  the  QAP  formulation  of  the  1-sum  are  both 
of  rank  n.  Hence  the  bounds  we  obtain  for  this  problem  by  this  approach  are  considerably 
more  involved. 

5.  Eigenvalue  bounds  for  the  2-sum  problem. 

5.1.  Orthogonal  bounds.  We  have  expressed  the  2-sum  problem  as  a  QAP  involving 
the  Laplacian  matrix  Q  and  the  matrix  B  —  where  p  is  an  n- vector  with  i-th  component 
equal  to  i.  The  trace  formulation  of  the  QAP  is 

imntr  QXBX^, 

where  X  ranges  over  the  set  of  permutation  matrices. 

A  technique  for  obtaining  lower  (upper)  bounds  for  a  QAP  is  to  relax  the  requirement 
that  the  minimum  (maximum)  be  attained  over  the  class  of  permutation  matrices.  Denote 
the  n-vector  of  all  ones  by  u.  A  matrix  X  of  order  n  is  a  permutation  matrix  if  and  only  if 
it  satisfies  the  following  three  constraints: 


(5.1) 

Xu 

=  u,  X^u  =  u: 

(5.2) 

X^X 

~  Inj 

(5.3) 

Xij 

>  0,  i,i  =  l,... 

The  first  of  these,  the  stochasticity  constraint,  expresses  the  fact  that  each  row  or  column  of 
a  permutation  matrix  has  a  single  nonzero  element  with  value  one;  the  second  states  that  a 
permutation  matrix  is  orthonormal;  and  the  third  that  its  elements  are  non-negative.  The 
simplest  bounds  for  a  QAP  are  obtained  when  we  relax  both  the  stochasticity  and  non¬ 
negativity  constraints,  and  insist  only  that  X  be  orthonormal.  The  following  result  is  from 
[10];  see  also  [11]. 

Theorem  5.1.  Let  the  eigenvalues  of  a  matrix  be  ordered 

Ai(-)<A2(-)---<A„(-). 

Then,  as  X  varies  over  the  set  of  permutation  matrices,  the  following  upper  and  lower  bounds 
hold: 

\i{Q)\n+i-i{B)  <  tr  QXBX^  <  k{Q)Xi{B).  □ 

i=l  1=1 

The  Laplacian  matrix  Q  has  Ai((5)  =  0;  also  Xi{B)  —  0,  for  i  =  1,  ...,  n  —  1,  and 
Xn{B)  —  ^ p=:  (l/6)n(n  -|-  l)(27z  -|- 1).  Hence  the  lower  bound  in  the  theorem  above  is  zero, 
and  the  upper  bound  is  {1  / 6) Xn{Q)n{n  4-  l)(2n  -|- 1). 

5.2.  Projection  bounds.  Stronger  bounds  can  be  obtained  by  a  projection  technique 
described  by  Hadley,  Rendl,  and  Wolkowicz  [17].  The  idea  here  is  to  satisfy  the  stochasticity 
constraints  in  addition  to  the  orthonormality  constraints,  and  relax  only  the  non-negativity 
constraints.  This  technique  involves  projecting  a  permutation  matrix  X  into  a  subspace 
orthogonal  to  the  stochasticity  constraints  (5.1)  by  means  of  an  eigen-projection. 
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From  now  on  we  normalize  u  :=  {\jy/n)u^  and  let  the  n  x  n  —  1  matrix  V  be  an 
orthonormal  basis  for  the  orthogonal  complement  of  u.  By  choice  of  V,  we  have  V^u  —  0, 
and  P  —  (^  u  y  ^  is  an  orthonormal  matrix  of  order  n. 

Observe  that 


P-XP  =  ( ^;  )  u  V) 

where  Y  =  V^XV. 

This  suggests  that  we  take 

X  = 

(5.4) 


(  iJXu  iJXV  \  /  1  0^  \ 
\  v'^Xu  v^xv  y  V  0  Y  ) 


uu^  +  VYV^. 


Note  that  with  this  choice,  the  stochasticity  constraints  Xu  —  u,  and  X'^u  —  u  are  satisfied. 
Furthermore,  if  X  is  an  orthonormal  matrix  of  order  n,  then 


P'^XP 


1  \ 
0  Y  ) 


is  orthonormal,  and  this  implies  that  Y  is  an  orthonormal  matrix  of  order  n  —  1.  Conversely, 
if  Y  is  orthonormal  of  order  n  —  1,  then  the  matrix  X  obtained  by  the  construction  above 
is  orthonormal  of  order  n.  The  non- negativity  constraint  X  >  0  becomes,  from  (5.4), 
VYV'^  >  —u  u^.  These  facts  will  enable  us  to  express  the  original  QAP  in  terms  of  a 
projected  QAP  in  the  matrix  of  variables  Y. 

To  obtain  the  projected  QAP,  we  substitute  the  representation  of  X  from  (5.4)  into  the 
objective  function  tr  QXBX^.  Since  Qu  =  0  by  the  construction  of  the  Laplacian,  terms  of 
the  form  Qu  u^  ■  ■  ■  vanish.  Further, 

tr  QVYV'^Bux/  =  tru^QVYV^u, 


where  we  use  the  identity  tr  MN  =  tr  NM  for  &n  n  x  k  matrix  M  and  a.  k  x  n  matrix 
N.  Again  this  term  is  zero  since  ujQ  =  0^.  Hence  the  only  nonzero  term  in  the  objective 
function  is 


tr  Q  VYV^  B  VY'^V^ 

=  triV^QV)Y  {V^BV)Y^ 

=  trQYBY^, 

where  M  =  VMV^. 

We  have  obtained  the  projected  QAP  in  terms  of  the  matrix  Y  of  order  n  —  1,  where 
the  constraint  that  X  be  a  permutation  matrix  now  imposes  the  constraints  that  Y  is  or¬ 
thonormal  and  that  VYV^  >  —uu^.  We  obtain  lower  and  upper  bounds  in  terms  of  the 
eigenvalues  of  the  matrices  Q  and  B  by  relaxing  the  non-negativity  constraint  again. 
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Theorem  5.2.  The  following  upper  and  lower  bounds  hold  for  the  2-sum  problem: 
{lll2)\2{Q){n  -  l)n{n  +  1)  <  <  {\l\2)\n{Q){n  -  l)n(n  +  1). 

Proof  If  we  apply  the  orthogonal  bounds  to  the  projected  QAP,  we  get 

I;  Xi{Q)K-i{B)  <  4(A)  <  i:  Ai(g)\(B). 

i-1  i=l 

The  vector  u  is  the  eigenvector  of  Q  corresponding  to  the  zero  eigenvalue,  and  hence  eigen¬ 
vectors  corresponding  to  higher  Laplacian  eigenvalues  are  orthogonal  to  it.  Thus  any  such 
eigenvector  Xj  can  be  expressed  as  Xj  =  Vr^.  Substituting  this  last  equation  into  the  eigen¬ 
value  equation  Qxj  =  Xj{Q)xj,  and  pre- multiplying  by  we  obtain  Qr^  =  Xj(Q)ry  Hence 
for  z  =  2,  . . .,  n,  we  have  Xi{Q)  —  A,_i((5).  Also,  A„_i(H)  =  •jfVV'^p,  and  all  other  eigen¬ 
values  are  zero.  Hence  it  remains  to  compute  the  largest  eigenvalue  of  B. 

From  the  representation  =  PP^  =  u  +  VV'^,  we  compute 

^VV^p  =  £  —  (p^  u)  p) 

=  (l/6)n(n  +  l)(2n  +  l)-(l/4)n(n  +  l)2  ^  (l/12)(n  -  l)n(n  +  1). 

We  get  the  result  by  substituting  these  eigenvalues  into  the  bounds  for  the  2-sum.  □ 

For  justifying  the  spectral  algorithm  for  minimizing  the  2-sum,  we  observe  that  the 
lower  bound  is  attained  by  the  matrix  Xq  =  -h  VRS^V^,  where  R  (S)  is  a  matrix 
of  eigenvectors  of  Q  (H),  and  the  eigenvectors  correspond  to  the  eigenvalues  of  Q  {B)  in 
increasing  (decreasing)  order. 

The  result  given  above  has  been  obtained  by  Juvan  and  Mohar  [22]  without  using  a  QAP 
formulation  of  the  2-sum.  We  have  included  this  proof  for  two  reasons:  First,  in  the  next 
subsection,  we  show  how  the  lower  bound  may  be  strengthened  by  diagonal  perturbations 
of  the  Laplacian.  Second,  in  the  following  section,  we  consider  the  problem  of  finding  a 
permutation  matrix  closest  to  the  orthogonal  matrix  attaining  the  lower  bound. 

5.3.  Diagonal  perturbations.  The  lower  bound  for  the  2-sum  can  be  further  im¬ 
proved  by  perturbing  the  Laplacian  matrix  Q  by  a  diagonal  matrix  Diag{d^,  where  d  is  an 
n- vector,  and  then  using  an  optimization  routine  to  maximize  the  smallest  eigenvalue  of  the 
perturbed  matrix. 

Choosing  the  elements  of  d  such  that  its  elements  sum  to  zero,  i.e.,  vTd  —  0,  simplifies 
the  bounds  we  obtain,  and  hence  we  make  this  assumption  in  this  subsection.  We  begin  by 
denoting  Q{d)  —  Q  +  Diag{d)^  and  expressing 

f{X)  =  trQ  XBX^  =  tr  Q{d)XBX^  -  tr  Diag{d)XBX^ . 

The  second  term  can  be  written  as  a  linear  assignment  problem  (LAP)  since  one  of  the 
matrices  involved  is  diagonal.  Let  the  permutation  vector  a  =  Xp^  and  let  d^  denote  the 
n-vector  formed  from  the  diagonal  elements  of  B. 

n 

tr  Diag{d)XBX^  =  ^  di6o,(i),o((j)  =  tr  ddg^ X^ . 

i=l 
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We  now  proceed,  as  in  the  previous  subsection,  to  obtain  projected  bounds  for  the 
quadratic  term,  and  thus  for  f{X).  Note  that  Q{d)u  =  {\ly/n)d  since  Qu  =  0;  and 
u^Q{d)u  =  0  since  the  elements  of  d  sum  to  zero.  We  shall  write  Bu  —  {\/y/n)r{B)  to 
denote  the  row-sum  of  the  elements  of  B. 

With  notation  as  in  the  previous  subsection,  we  substitute  X  —  ui^  +  VYV^  in  the 
quadratic  term  in  f{X).  The  first  term  tr  Q{d)uvjBui^  =  tr  u^Q{d)uu'^ Bu  =  0.  The 
second  and  third  terms  are  equal,  and  their  sum  can  be  transformed  as  follows: 

2  tr  Q{d)VYV^Bu  =  2  tr  u^Q{d)VYV^Bu 
=  {2ln)tr  £VYV'^r{B)  =  {2ln)  tr  V^r{B)  d^VY 
=  {2ln)trY^V^dr{BfV  =  {2/n)  tr  dr{BfVY^V^ . 

Note  that  this  term  is  linear  in  the  projected  variables  Y,  and  we  shall  find  it  convenient  to 
express  it  in  terms  of  X  by  the  substitution  X^  —  n  =  VY'^V^.  Thus 

{2/n)trdr{BfVY^V^  =  i2/n)  tr  dr{Bf{X^  -  uu^)  =  {2/n)  tr  dr{Bf  X^ , 

since  the  second  term  is  equal  to  tr  dr{B)^u^  which  is  zero  by  the  choice  of  d. 

Finally,  the  fourth  term  becomes  tr  Q{d)Y BY^ ,  where  Q{d)  =  V'^Q{d)V,  and  as  before 
B  =  V'^BV. 

Putting  it  all  together,  we  obtain 

fiX)  =  tr  Qid)YBY^  +  tr  {{2/n)V^dr{Bfx'^  -dde^X^)  . 

Observe  that  the  first  term  is  quadratic  in  the  projected  variables  Y,  and  the  remaining  terms 
are  linear  in  the  original  variables  X.  Our  lower  bound  for  the  2-sum  shall  be  obtained  by 
minimizing  the  quadratic  and  linear  terms  separately. 

We  can  simplify  the  linear  assignment  problem  by  noting  that  B  =  pp^.  Thus  rB,i  = 
—  (l/2)n(n  +  l)z,  and  hence  {2/n)r{B)  =  (n  +  l)p.  Further,  dg  =  sq[p),  the  vector 
with  zth  component  equal  to  B .  Hence  the  final  expression  for  the  linear  assignment  problem 
is 


tr  d  (^{n  -{-  l)p^  -  sq{p)^^  X^ . 

Let  L  denote  the  minimum  value  of  this  problem  (over  the  set  of  permutation  matrices  X, 
for  a  given  d),  which  can  be  computed  as  the  solution  of  a  transportation  problem. 

The  eigenvalues  of  B  can  be  computed  as  in  the  previous  subsection.  We  may  choose  d 
to  maximize  the  smallest  eigenvalue  of  the  matrix  Q{d).  Thus  this  discussion  leads  to  the 
following  result. 

Theorem  5.3.  The  minimum  2-sum  of  a  symmetric  matrix  A  can  be  bounded  as 

>  max(l/12)Ai((5(d))(n  -  l)n(n  +  1)  +  X, 
d 

where  the  components  of  the  vector  d  sum  to  zero.  □ 
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6.  Computing  an  approximate  solution  from  the  lower  bound.  Consider  the 
problem  of  finding  a  permutation  matrix  Z  “closest”  to  an  orthogonal  matrix  Xq  that  attains 
the  lower  bound  in  Theorem  5.2.  We  show  in  this  section  that  sorting  the  second  Laplacian 
eigenvector  components  in  non-increasing  (also  non-decreasing)  order  yields  a  permutation 
matrix  that  solves  a  linear  approximation  to  the  problem.  This  justifies  the  spectral  approach 
for  minimizing  the  2- sum. 

Since  Xq  attains  the  lower  bound  in  Theorem  5.2,  Xq  =  u'lF  -\-  VRS'^V'^,  where  R  {S) 
is  a  matrix  of  eigenvectors  of  Q  {B)  corresponding  to  the  eigenvalues  of  Q  (B)  in  increasing 
(decreasing)  order.  We  begin  with  a  preliminary  discussion  of  some  properties  of  the  matrix 
Wo  and  the  eigenvectors  of  Q.  For  j  =  1,  . . .,  n  —  1,  let  the  jth  column  of  R  be  denoted 
by  Tj,  and  similarly  let  Sj  denote  the  jth  column  of  S.  Then  Sj  =  V^p,  and  for  j  =  2,  . . ., 
n  —  1,  the  vector  Sj  is  orthogonal  to  i.e., 

(6.1)  s/l/^p  =  0. 

Recall  from  the  previous  section  that  a  second  Laplacian  eigenvector  X2  —  ^T.\- 

Now  we  can  formulate  the  “closest”  permutation  matrix  problem  more  precisely.  The 
minimum  2-sum  problem  may  be  written  as 

min  \\{Q  +  aIfl‘^Zp\\2  ■ 

Zi 

We  have  chosen  a  positive  shift  a  to  make  the  shifted  matrix  positive  definite  and  hence 
to  obtain  a  weighted  norm  by  making  the  square  root  nonsingular.  It  can  be  verified  that 
the  shift  has  no  effect  on  the  minimizer  since  it  adds  only  a  constant  term  to  the  objective 
function. 

We  substitute  Z  —  Xo  +  {Z  —  Wo)  and  expand  the  2-sum  about  Wq  to  obtain 
IKO  +  = 

(6.2)  11(0  +  a/)>/W„£||2"  +  2ir  ^(Z  -  Xaf  (Q  +  aI)X„p  +  ||(0  +  alfl^iZ  - 

The  first  term  on  the  right-hand-side  is  a  constant  since  W  is  a  given  orthogonal  matrix; 
the  third  term  is  a  quadratic  in  the  difference  {Z  —  Wo)  and  hence  we  jieglect  it  to  obtain  a 
linear  approximation.  It  follows  that  we  can  choose  a  permutation  matrix  Z  close  to  Wq  to 
approximately  minimize  the  2-sum  by  solving 

(6.3)  mintr  p^Z^((5  +  o;/)Wo£  =  ram  tr  [Q aI)XQB Z^ . 

7j  Tj 

Substituting  for  Wo  from  (5.4)  in  this  linear  assignment  problem  and  noting  that  Qu  =  0, 
we  find 

min  tr  (Q  -F  aI)XoBZ^  =  min  tr  {Q  +  al)  (u  +  VRS^V^)BZ^ 

z  z 

(6.4)  =  mm  tr  QVRS^V^BZ^  -b  atruu^BZ^  +  atr  VRS^V^BZ^. 
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The  second  term  on  the  right-hand-side  is  a  constant  since 


truj/BZ'^ 

=  p)  tr  u 
=  (m^  p)  tr  Z^u 
=  p)  tr  u  jF 
=  (n^p)T 

Here  we  have  substituted  Z’^u  —  u  from  (5.1).  We  proceed  to  simplify  the  first  term  in  (6.4), 
which  is 


tr  QVRS'^BZ'^  =  tr  QV  f  g  r,  «/)  V'^pfZ^. 

From  (6.1)  we  find  that  s^V'^p  =  0,  for  ^’  =  2,  . . .,  n  —  1,  and  hence  only  the  first  term  in 
the  sum  survives.  Noting  that  s^  —  V^p,  and  V"ri  =  X2,  this  term  becomes 

tr  Qx2  p^VV'^pp^Z'^  =  >^2iQ){p^VV'^p)  tr  X2  ^ Z'^ . 

The  third  term  in  (6.4)  can  be  simplified  in  like  manner,  and  hence  ignoring  the  constant 
second  term,  this  equation  becomes 

('^2(Q)  +  Oi)  i^VV^p)  min  tr  X2  p^Z^. 

Hence  we  are  required  to  choose  a  permutation  matrix  Z  to  minimize  tr  X2  p'^ Z^  = 
tr  Z^ X2  p^.  The  solution  to  this  problem  is  to  choose  Z  to  correspond  to  a  permutation 
of  the  components  of  ^2  in  non-increasing  order,  since  the  components  of  the  vector  p  are 
in  increasing  order.  Note  that  —X2  is  also  an  eigenvector  of  the  Laplacian  matrix,  and 
since  the  positive  or  negative  signs  of  the  components  are  chosen  arbitrarily,  sorting  the 
eigenvector  components  into  non-decreasing  order  also  gives  a  permutation  matrix  Z  closest 
to  the  orthogonal  matrix  Xq. 

Similar  techniques  can  be  used  to  show  that  if  one  is  interested  in  maximizing  the  2-sum, 
then  a  closest  permutation  matrix  to  the  orthogonal  matrix  that  attains  the  upper  bound 
in  Theorem  5.2  is  approximated  by  sorting  the  components  of  the  Laplacian  eigenvector 
(corresponding  to  the  largest  eigenvalue  A„((5))  in  non-decreasing  (non-increasing)  order. 

7.  Implications.  In  this  section  we  consider  the  implications  of  the  spectral  lower 
bounds  that  we  have  obtained.  We  denote  the  eigenvalues  \2{Q)  by  A2  and  \n{Q)  by  A„ 
for  the  sake  of  brevity  in  this  section.  We  make  use  of  the  following  result  whose  proof 
may  be  found  in  [33].  Let  {G„}  represent  a  family  of  graphs  with  parameter  n  representing 
the  number  of  vertices,  and  let  G  be  a  graph  on  n  vertices  from  this  family.  Then  G  has 
an  /(n)-separator  S  if  |5|  =  ©(/(«)),  and  S  separates  G  \  S  into  two  parts  A,  B  with 
hn  <  |y4|,  |.B|  <  (1  —  /i)n,  where  0  <  h  <  1.  (Thus  A  and  B  have  0(n)  vertices.)  For 
instance,  planar  graphs  have  -y/n-separators. 
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Theorem  7.1.  If  the  family  {Gn}  has  f(n)-separators,  then 

\2  <  KA{n)^, 

n 

where  K  is  a  constant  independent  of  n,  and  A(n)  is  the  maximum  degree  of  any  vertex  in 

G.  □ 

Results  of  this  nature,  involving  separator  properties  of  graphs,  provide  the  tightest 
known  upper  bounds  for  the  second  eigenvalue.  Other  results  [8]  that  bound  A2  in  terms  of 
the  edge  or  vertex  connectivity  lead  to  an  upper  bound  equal  to  at  least  one. 

We  can  use  information  about  the  asymptotic  behavior  of  the  second  eigenvalues  together 
with  the  bounds  we  have  obtained  to  predict  the  behavior  of  envelope  parameters.  For  the 
envelope  size,  we  make  use  of  Theorem  3.2;  for  the  estimate  of  the  envelope  work  Wbound{A), 
we  will  need  to  combine  Theorem  5.2  and  Theorem  2.1  to  obtain  the  following  result. 

Theorem  7.2.  The  estimate  of  the  minimum  work  in  an  envelope  factorization  scheme 
of  a  symmetric  positive  definite  matrix  A  is  bounded  by 

WboundminiA)  >  ^  X2{Q)  (u  -  l)n(n  +  1).  □ 

(Note  that  we  have  chosen  to  use  the  simpler  bounds  in  Theorem  5.2  rather  than  the 
bound  from  Theorem  5.3.) 

The  bound  on  envelope  size  is  tight  for  rather  dense  graphs  and  matrices.  For  instance, 
the  full  matrix  (the  complete  graph)  has  A2  =  A  +  1  =  n,  and  hence  Esizem.in{A.)  =  0(n^). 
Similarly  we  find  that  the  bound  on  the  envelope  work  Wboundmini-^)  =  Moreover, 

the  predicted  lower  bound  is  within  a  factor  of  three  of  the  envelope  size.  These  bounds  are 
also  asymptotically  tight  for  random  graphs  where  each  possible  edge  is  present  in  the  graph 
with  a  fixed  probability  p,  since  the  second  Laplacian  eigenvalue  satisfies  [21] 

X2  =  pn-  0([p(l  -p)n  logn]^/^). 

We  are  interested  in  understanding  the  implications  of  these  bounds  for  finite  element 
meshes  in  two  and  three  dimensions,  with  the  maximum  degree  A  bounded  independent 
of  n.  For  planar  meshes,  since  /(n)  =  (!7(n^/^),  we  have  that  A2  =  Note  that 

‘well-shaped’  two-dimensional  finite  element  meshes  have  A2  =  0(/i^)  ^  0(n“^),  where  h  is 
the  (uniform)  distance  between  two  successive  points  of  the  mesh  along  the  x-  or  j/-axis.  For 
three-dimensional  meshes  that  possess  separators,  we  have  A2  =  Again  for 

‘good’  three-dimensional  meshes,  we  have  A2  =  0(/i^)  =  0(n“^/^). 

Plugging  in  the  bounds  above  for  A2  into  the  envelope  size  and  envelope  work  bound, 
we  find  that  for  planar  meshes  that  have  A2  =  0(n“^'^^),  the  envelope  size  Esizemin{A)  = 
fl(n^/^),  and  the  work  bound  is  Wboundmin{A)  =  fi(n®/^).  For  planar  meshes  with  A2  = 
0(n~^),  we  obtain  the  lower  bounds  Esizemin{A)  =  0(n),  and  WboundminiX^)  =  fl(n^).  For 
three-dimensional  meshes  that  have  A2  =  0(n"^/^),  we  have  Esizemin{A)  =  fl(n®/^)  and 
Wboundmin{A)  =  n(n^/^).  Again  for  ‘good’  three-dimensional  meshes  with  A2  =  0(n~^/^), 
we  have  EsizCmini-d^)  =  and  Wboundmin{A)  =  fi(n^/^). 

It  is  easy  to  establish  that  A„  <  2A  for  all  graphs,  by  considering  the  eigenvalue  equa¬ 
tion  Qx^  =  A„x„.  Since  we  are  considering  bounded- degree  finite  element  meshes,  from 
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Theorem  3.2,  Theorem  5.2,  and  Theorem  2.1,  the  upper  bound  on  the  envelope  size  is 
Esize{A)  —  0{n'^)  and  that  on  the  work  bound  is  Wbound{A)  =  0{n^).  Notice  that  the 
gap  between  the  lower  bound  and  the  upper  bound  on  these  parameters  is  quite  large. 

How  good  are  the  lower  bounds  on  the  envelope  size  and  estimate  of  the  envelope  work? 
To  answer  this  question,  we  prove  upper  bounds  for  these  parameters  from  a  “modified 
nested  dissection”  ordering  of  graphs  wdth  good  separators.  At  each  level,  if  a  separator  S 
separates  the  current  graph  G  into  two  parts  A,  B,  then  in  this  ordering  vertices  in  A  are 
numbered  first,  then  the  vertices  in  S,  and  finally  the  vertices  in  B.  (See  the  ordering  in 
Figure  2.1,  where  S  corresponds  to  the  set  of  vertices  in  the  middle  column.)  For  simplicity 
of  exposition,  initially  consider  the  class  of  planar  graphs. 

Theorem  7.3.  The  minimum  envelope  size  of  a  symmetric  matrix  of  order  n,  M ,  whose 
adjacency  graph  G  is  a  bounded-degree  planar  graph,  satisfies  Esize,nin{AI)  =  0{n^-^). 

Proof  Let  a  separator  S  separate  G  into  parts  A  and  B,  with  |A|  —  an  >  \B\  =  /?n,  and 
151  =  c^/n,  where  a-\-  ^  =  1.  From  the  planar  separator  theorem  [26],  a  is  bounded  between 
1/2  and  2/3.  Consider  the  modified  nested  dissection  ordering  that  numbers  vertices  in  A 
first,  then  vertices  in  S,  and  finally  vertices  in  B.  The  contribution  to  the  envelope  size 
made  by  vertices  in  S  is  at  most  Cy/nan,  since  each  of  them  has  row-width  at  most  |A|. 
Also,  there  are  at  most  Ac^/n  vertices  in  B  adjacent  to  vertices  in  5;  each  of  these  vertices 
has  row-width  at  most  fdn  since  they  may  be  numbered  last  in  B,  and  hence  such  vertices 
contribute  at  most  Cy/nA/Sn  to  the  envelope  size. 

The  reasoning  above  leads  to  the  recurrence 

E{n)  <  E{an)  +  E{(3n)  -f  c(a  +  /?A) n^'^. 

Let  a  =  a^-^  Then  it  is  easy  to  show  that  a  <  1  since  q  +  /3  =  1.  The  solution  to  the 

recurrence  is 

E{n)  <  c(a  +  ^A)  — - — n^'^. 

1  —  G 

This  completes  the  proof.  □ 

This  technique  can  be  used  to  prove  that  Wbound{A)  —  for  matrices  whose 

adjacency  graphs  are  bounded-degree  planar  graphs.  For  matrices  whose  adjacency  graphs 
(i)  have  bounded-degree,  (ii)  are  embeddable  in  three  dimensions,  and  (iii)  possess 
separators,  the  upper  bounds  are  Esize{A)  =  and  Wbound{A)  =  0{n^^^).  These 

results  show  that  the  lower  bounds  on  the  minimum  envelope  size  (and  the  estimate  of  the 
work-bound)  are  asymptotically  tight  for  two-dimensional  meshes  with  A2(G)  =  0(n"^/^), 
and  for  three-dimensional  meshes  with  A2(G)  =  0(n"^/^).  Upper  bounds  as  in  the  previous 
theorem  can  also  be  established  for  the  class  of  bounded-degree  overlap  graphs  embedded 
in  d-dimensional  space  [28];  these  graphs  have  0{{d  -  l)/d)-separators  that  separate  them 
into  two  parts  with  at  most  ((d  -t-  l)/(d  -f  2))n  vertices  in  the  larger  part. 

The  broader  implication  of  the  result  above  is  that  if  a  problem  possesses  good  separators, 
then  it  has  small  envelope  size  and  work.  In  practice,  there  exist  orderings  with  smaller 
envelope  parameters  than  the  modified  nested  dissection  ordering. 

We  conclude  that  for  problems  in  dimensions  greater  than  two,  asymptotically  envelope 
factorization  schemes  are  not  competitive  with  iterative  methods  or  direct  methods  that 
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18 

48 

2.00 

969 

978 

0.9 

66 

192 

6.25e-l 

1.50e-|-4 

1.54e+4 

2.6 

258 

768 

1.65e-l 

2.36e+5 

2.53e+5 

6.9 

1,026 

3,072 

4.17e-2 

3.75e+6 

4.05e+6 

7.4 

4,098 

12,270 

1.05e-2 

6.00e-k7 

6.44e-h7 

7.3 

16,386 

49,152 

2.60e-3 

0.953ed-9 

1.03e+9 

9.1 

Table  8.1 


2-sums  from  the  spectral  reordering  algorithm  and  lower 


bounds  for  triangulations  of  the  sphere. 


store  and  operate  on  only  the  nonzero  elements.  (For  the  latter,  the  storage  requirements 
are  O(nlogn)  and  work  is  for  a  two-dimensional  problem;  in  three  dimensions, 

these  are  ©(ra^/^logn)  and  0(n^).)  But  when  a  two-dimensional  mesh  possesses  a  small 
second  Laplacian  eigenvalue,  envelope  methods  may  be  expected  to  work  well.  Similar 
conclusions  should  hold  for  three-dimensional  problems  when  the  number  of  mesh-points 
along  the  third  dimension  is  small  relative  to  the  number  in  the  other  two  dimensions,  and 
for  two-dimensional  surfaces  embedded  in  three-dimensional  space. 

8.  Computational  results.  We  present  computational  results  to  verify  how  well  the 
spectral  ordering  reduces  the  2-sum.  We  report  results  on  two  sets  of  problems. 

The  first  set  of  problems,  shown  in  Table  8.1,  is  obtained  from  John  Richardson’s  (Think¬ 
ing  Machines  Corporation)  program  for  triangulating  the  sphere.  The  spectral  lower  bounds 
reported  are  from  Theorem  5.2.  The  results  show  that  the  spectral  reordering  algorithm 
computes  values  within  a  few  percent  of  the  optimal  2-sum,  since  the  gap  between  the  2-sum 
lower  bounds  and  the  spectral  2-sums  is  within  that  range. 

Table  8.2  contains  the  second  set  of  problems,  taken  from  the  Boeing-Harwell  and  Nasa 
collections.  Here  the  bounds  are  weaker  than  the  bounds  in  Table  8.1.  These  problems 
have  two  features  that  distinguish  then  from  the  sphere  problems.  Many  of  them  have 
less  regular  degree  distributions — e.g.,  NASA1846  has  maximum  degree  41  and  minimum 
degree  5,  though  NACA  is  an  exception.  They  also  represent  more'complex  geometries. 
Nevertheless,  these  results  imply  that  the  spectral  2-sum  is  within  a  factor  of  two  of  the 
optimal  value  for  these  problems. 

The  gap  between  the  computed  2-sums  and  the  lower  bounds  could  be  further  reduced 
in  two  ways.  First,  a  local  reordering  algorithm  applied  to  the  ordering  computed  by  the 
spectral  algorithm  might  potentially  decrease  the  2-sum.  Second,  the  lower  bounds  could  be 
improved  by  incorporating  diagonal  perturbations  to  the  Laplacian.  We  will  consider  both 
these  issues  in  future  work. 

9.  Conclusions.  The  lower  bounds  on  the  2-sums  show  that  the  spectral  reordering 
algorithm  can  yield  nearly  optimal  values.  To  the  best  of  our  knowledge,  these  are  the  first 
results  providing  reasonable  bounds  on  the  quality  of  the  orderings  generated  by  a  reordering 
algorithm  for  minimizing  envelope-related  parameters.  Earlier  work  had  not  addressed  the 
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NASA1824 
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1.37e+8 

1.74e+8 

21 
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2,146 

35,052 

1.35e-l 

l.lle+8 

1.32e+8 

16 

NACA 

4,224 

12,416 

3.57e-3 

2.24e-b7 

2.70e-f7 

17 

BARTH4 

6,019 

17,473 

1.756-3 

3.19e+7 

5.41e+7 

41 

BARTH 

6,691 

19,478 

2.62e-3 

6.55e+7 

1.39e+8 

53 

Table  8.2 

2-sums  from  ihe  spectral  reordering  algorithm  and  lower  bounds  for  some  problems  from  the  Boeing- 
Harwell  and  Nasa  collections. 


issue  of  the  quality  of  the  orderings  generated  by  the  algorithms.  The  lower  bounds  for  the 
envelope  size  and  envelope  work  are  not  useful  in  assessing  the  quality  of  envelope-reducing 
orderings,  but  are  asymptotically  tight  for  certain  classes  of  2-  and  3-dimensional  meshes.  It 
remains  to  be  seen  if  lower  bounds  obtained  via  the  QAP  formulation  of  the  1-sum  problem 
lead  to  better  results  for  the  envelope  size.  Designing  local  reordering  algorithms  to  improve 
the  quality  of  the  spectral  orderings  is  another  priority. 
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Appendix 


A.  Lower  bounds  on  the  minimum  p-sum.  We  prove  two  lower  bounds  on  the 
minimum  p-sums.  We  make  use  of  Lemma  3.1  in  proving  the  first  result.  In  the  following 
Bm{x)  is  the  mth  Bernoulli  polynomial,  and  Bm  is  the  mth  Bernoulli  number. 

Theorem  A.l.  For  1  <  p  <  oo,  the  minimum  p-sum  of  a  graph  G  on  n  vertices 
satisfies 

>  {Bp+i{s  +  1)  -  .Bp+i) , 

p  +  1 

where  s  =  {X2f4:A)n. 

Proof.  Consider  any  ordering  a  of  the  vertices  of  G.  Partition  the  vertices  into  two  sets: 
A  consisting  of  the  lowest-numbered  n/2  vertices,  and  B  consisting  of  the  highest-numbered 
n/2  vertices.  By  Lemma  3.1  the  number  of  edges  joining  A  and  B,  |^(A,B)|,  is 

|i(A,B)|>h(„/2f. 

n 

Hence  at  least  s  =  |(5(A,  5)|/A  vertices  in  B  are  adjacent  to  vertices  in  A.  Each  vertex 
in  this  subset  of  B  has  the  least  row- width  when  it  is  adjacent  to  the  highest-numbered 
vertex  in  A  and  to  no  other  vertices  in  A.  Hence  these  s  vertices  make  a  contribution  of  at 
least  P  -f  . . .  +  sP  to  the  p-sum,  and  this  sum  can  be  expressed  in  terms  of  the  Bernoulli 
polynomials  as  stated.  □ 

From  an  expansion  of  the  Bernoulli  polynomial,  we  find  that  asymptotically 

^  +  0((A27A’’)n'). 

We  proceed  to  obtain  another  lower  bound  on  the  minimum  p-sum. 

The  next  result  makes  use  of  the  following  Lemma  A. 2  recently  proved  by  Helmberg  et 
al.  [18].  Define  the  following  symmetric  function  of  the  two  positive  integers  mi,  m2  (with 
mi  -f  m2  <  n)  and  parameters  A2,  A„: 


(A.l)  /(mi,m2)  = 

— ^  [ymim2  +  \J{n-  mi)(n  -  m2)^  A2  -f  {yfnffrvi  -  \J{n-  mfijin  -  m2)^  A,. 

Lemma  A. 2.  Let  Si,  S2  be  two  disjoint  subsets  of  the  vertices  of  a  graph  G  onn  vertices, 
with  l^il  =  Si,  for  i  =  2.  Then  the  number  of  edges  joining  Si  and  S2,  1^(51,52)1,  satisfies 


|^(5i,52)|  >  /(5i,S2).  □ 
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Theorem  A. 3.  For  1  <  p  <  oo,  the  minimum  p-sum  of  a  graph  G  satisfies 


73  \  ^ 

""p  ^  2P+1A  (A„  +  A2)p+2 


(2A„  4-  A2)(A„  +  2A2)w^^^. 


Proof.  Consider  any  ordering  a  of  the  vertices  of  G,  and  consider  a  tripartition  A,  B, 
C:  We  choose  A  to  consist  of  the  lowest-numbered  a  =  (n  —  b)/2  vertices,  C  to  consist  of 
the  highest-numbered  (n  —  b)/2  vertices,  and  B  to  contain  the  remaining  b  vertices  in  the 
‘middle’.  Here  6,  the  size  of  jB,  is  a  parameter  that  will  be  determined  later  to  obtain  a  large 
lower  bound. 

From  Lemma  A. 2,  |6(A,  (7)1,  the  number  of  edges  joining  A  and  (7,  is  at  least  /(a,  a), 
where  the  symmetric  function  /(•,.)  is  defined  in  (A.l).  Hence  there  are  at  least  sc  = 
f{a,a)IA  vertices  in  C  adjacent  to  vertices  in  A.  Each  of  these  vertices  has  row- width  at 
least  b. 

Initially,  consider  the  contribution  to  the  envelope  size  Esize[G)  made  by  these  vertices 
to  obtain  a  suitable  value  for  b. 

(A. 2)  Esize(G)  > 


We  choose  b  to  maximize  the  lower  bound  on  the  envelope  size.  Differentiating  the  cubic 
polynomial  in  (A. 2)  with  respect  to  b  and  simplifying,  we  obtain  the  quadratic  equation 

2  2  A2  +  A„  1 A2  2  n 

b  -- — - nb+-—n  =0. 

O  Aji  O  An 


A 

(n-b) 


in 


'  n  —  b  n  +  b\  . 
^r-  +  A2  + 


'n  —  b  n  b'' 


1 

4A 


6(n  -  6)  (A2  -  (6/n)A„) 


From  the  quadratic  we  find  that  the  maximizer  is,  to  first  order,  =  (1/2)(A2/(A„  -|-  A2))n. 

Now  we  consider  the  contribution  to  the  p-sum  made  by  the  sc  vertices  in  C  adjacent 
to  vertices  in  A.  Each  of  these  vertices  contributes  at  least  If  to  the  p-sum,  and  thus  a  lower 
bound  on  the  minimum  p-sum  is 

<(<?)>  J^(n-*>)(A2-(6/n)A„)6'. 

It  is  not  easy  to  find  a  maximizer  of  the  right-hand-side  in  the  bound  above  on  the  p-sum 
since  the  polynomial  in  b  is  of  degree  p  +  2.  Hence  we  choose  b  equal  to  the  maximizer  of  the 
envelope  size.  We  obtain  the  bound  stated  in  theorem  by  substituting  6  =  6m  in  the  bound 
above.  □ 

Juvan  and  Mohar  [22]  have  proved  upper  bounds  for  the  p-sums.  The  techniques  in  this 
Appendix  can  be  used  to  compute  bounds  on  Esize{A)  and  Wbound{A),  but  the  results  are 
weaker  than  those  obtained  in  Section  3. 
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